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ABSTRACT 

The extended CASA Line Analysis Software Suite (XCLASS) is a toolbox for the Common Astronomy Software Applications 
package (CASA) containing new functions for modeling interferometric and single dish data. Among the tools is the my XCLASS 
program which calculates synthetic spectra by solving the radiative transfer equation for an isothermal object in one dimension, 
whereas the finite source size and dust attenuation are considered as well. Molecular data required by the myXCLASS program are 
taken from an embedded SQLiteS database containing entries from the Cologne Database for Molecular Spectroscopy CDMS) and 
JPL using the Virtual Atomic and Molecular Data Center (VAMDC) portal. Additionally, the toolbox provides an interface for the 
model optimizer package Modeling and Analysis Generic Interface for eXtemal numerical codes (MAGIX), which helps to find the 
best description of observational data using myXCLASS (or another external model program), that is, finding the parameter set that 
most closely reproduces the data. 
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1. Introduction 

The CASA package (McMullin et al. 2007) provides a powerful 
tool for data post-processing for ALMA and JVLA, but contains 
only rudimentary functions for modeling the data. But modeling 
astronomical data is essential to derive physical parameters such 
as column densities and rotational temperatures, as well as in¬ 
formation on the location and the kinematics of the emitting gas 
component for each observed molecular species. 

The toolbox described in this paper offers not only the 
possibility to model data by using the myXCLASS program 
(Schilke et al. 2001; Comito et al. 2005; Zernickel et al. 2012; 
Crockett et al. 2014a; Crockett et al. 2014b; Neill et al. 2014), 
which computes synthetic spectra by solving the detection equa¬ 
tion similar to software packages like Weeds (Maret et al. 201 1) 
and CASSIS* *, but also makes use of the optimization package 
MAGIX (Molleref a/. 2013). MAGIX provides a framework of 
an easy interface between existing codes and an iterating en¬ 
gine that attempts to minimize deviations of the model results 
from available observational data, constraining the values of the 
model parameters and providing corresponding error estimates. 
In addition to the myXCLASS program many external model pro¬ 
grams^ such as RADEX (Van der Tak et al. 2007), RADMC- 
3D* or LIME (Brinch & Hogerheijde 2010) can be plugged into 
MAGIX to explore their parameter space and find the set of pa¬ 
rameter values that best fits observational data. Most of the opti- 


* http://WWW.astro.uni-koeln.de/proj ects/schilke/ 
myXCLASSInterface 
' http://cassis.irap.omp.eu/ 

^ Here, the phrase “external model program” means the external pro¬ 
gram that calculates the model function depending on several input pa¬ 
rameters. 

^ http: //WWW. ita.uni-heidelberg.de/~dullemond/ 

software/radmc-3d/ 


mization algorithms included in the MAGIX package are available 
as an OpenMP"* and an MPI^ version to speed up the computa¬ 
tion. Eurthermore, the toolbox includes two new functions which 
provide a simplified interface for MAGIX in conjunction with 
myXCLASS to model single spectra and complete data cubes, re¬ 
spectively. Once installed, it is seamlessly integrated into CASA, 
so the tasklist command within CASA lists all the new func¬ 
tions. In addition, the help command followed by the name of 
the function provides a short description of the corresponding 
function with a short overview of all input parameters. In addi¬ 
tion, the XCLASS interface can be used without CASA, where¬ 
fore the following python packages has to be installed: numpy, 
scipy, pyfits, matplotlib and sqliteB. 

The Splatalogue database® is accessible from within the 
viewer in CASA. We have opted to use our own SQLite3 
database based on the VAMDC^ implementation, mostly be¬ 
cause this allows the user to augment the data by providing the 
partition function between 1.072 and 1000 K at 110 different 
temperatures in contrast to seven temperatures in the Splatalogue 
database, which is derived from the traditional JPL* and CDMS® 
catalogs. It is hoped that the content of Splatalogue and the 
VAMDC database will at some point be homogenized. 

The XCLASS interface for CASA is mostly written in 
python, whereas the myXCLASS program and most of the algo¬ 
rithms included in the MAGIX package are written in Lortran 90. 
All required Python modules are already included in the CASA 
package, so no further software package except the GNU com- 


http://openmp.org/wp/ 

^ http://www.open-mpi.org/ 

® http://www.splatalogue.net 
^ http://www.vamdc.eu 

* http://spec.jpl.nasa.gov 
^ http://www.cdms.de 
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piler'° (gfortran and gcc) with OpenMP/OpenMPI extension 
has to be installed. The XCLASS interface for CASA is avail¬ 
able for Linux and Mac OS 10.11 (64 bit). 

In the next sections we will describe the myXCLASS pro¬ 
gram and its algorithms, assumptions, and shortcomings in de¬ 
tail. Afterwards, we describe the SQLite3 database used for 
the myXCLASS program and the related functions provided by 
the XCLASS interface. In the following we briefly present the 
MAGIX package with a short overview of the included algo¬ 
rithms. Finally, we give a detailed description of the new func¬ 
tions myXCLASSFit and myXCLASSMapFit used for modeling 
data with myXCLASS and MGIX. 


2. myXCLASS 

The XCLASS interface for CASA contains the myXCLASS pro¬ 
gram, originally developed by P. Schilke (SchilkeefaL 2001; 
Comito et al. 2005; Zernickel et al. 2012; Crockett et al. 2014a; 
Crockett et al. 2014b; Neill et al. 2014) based on the GILDAS'^ 
package CLASS, which models a spectrum by solving the ra¬ 
diative transfer equation for an isothermal object in one dimen¬ 
sion, the detection equation (Stabler & Palla 2005). Here, LTE 
is assumed, that is, the source function is given by the Planck 
function of an excitation temperature, which does not need to 
be the physical temperature, but is constant for all transitions. 
The myXCLASS function is designed to describe line-rich sources 
which are often dense, so that LTE is a reasonable approxima¬ 
tion. Also, a non-LTE description requires collision rates which 
are available only for a few molecules. 

The myXCLASS function is able to model a spectrum with an 
arbitrary number of molecules where the contribution of each 
molecule is described by an arbitrary number of components. 
The 1-d assumption imposes a very simplistic geometrical struc¬ 
ture. We recognize two classes of components. 

One, the core objects (in earlier implementations called 
emission component), consists of an ensemble of objects cen¬ 
tered at the middle of the beam. These could be identified with 
clumps, hot dense cores etc. which overlaps but do not inter¬ 
act either because they do not overlap in physical or in veloc¬ 
ity space. For computational convenience, they are assumed to 
be centered in the beam, as shown in Fig. 1 . It is also assumed 
that the dust emission emanates (partly) from these components. 
Their intensities are added, weighted with the beam filling fac¬ 
tor, see Eq. (2). 

The second class, foreground objects (in earlier implemen¬ 
tations called absorption components), are assumed to be in lay¬ 
ers in front of the core components. In the current 1-d imple¬ 
mentation, they would have a beam averaged intensity of the 
core sources as background, and would fill the whole beam. 
Examples for such structures would be source envelopes in front 
of dense cores, or absorption components along the line-of-sight. 

As shown in Eig. 1, we assume that core components do not 
interact with each other radiatively, that is, one core layer is not 
influenced by the others. But the core layers may overlap to offer 
the possibility to model sources consisting of several molecules 
and compounds. The solution of the radiative transfer equation 


https://gcc.gnu.org/ 

" http://irain.fr/IRAHFR/GILDAS/ 



beam 


Fig. 1. Sketch of a distribution of core layers within the Gaussian 
beam of the telescope (black ring). Here, we assume three dif¬ 
ferent core components 1, 2, and 3, centered at the middle of 
the beam with different source sizes, excitation temperatures, 
velocity offsets etc. indicated by different colors. Additionally, 
we assume that all core components have the same distance to 
the telescope, i.e., all core layers are located within a plane per¬ 
pendicular to the line of sight. Furthermore, we assume that this 
plane is located in front of a background Layer 4 with homoge¬ 
neous intensity /^™^(v) over the whole beam. Core components 
do not interact with each other radiatively. 


for core layers is*^, 

+ - ./cmb) > (1) 

where the sums go over the indices m for molecule, and c for 
(core) component, respectively. In the following we will briefly 
describe each term in Eq. (1). 




The beam filling (dilution) factor of molecule m and 

component c in Eq. (1) for a source with a Gaussian brightness 
profile, see below, and a Gaussian beam is given by'^ 




(^.c)2 


6»,(v)2 H-(0'"’O2’ 


( 2 ) 


where and 6, represents the source and telescope beam full 
width half maximum (EWHM) sizes, respectively. The sources 
beam EWHM sizes for the different components are defined 
by the user in the molfit file, described in Sect. 2.1 . Additionally, 
we assume for single dish observations, that the telescope beam 
EWHM size is related to the diameter of the telescope by the 
diffraction limit 

fl,(v) = (L22-^)-^ = (l.22-(3) 


*2 A derivation of the expression can be found in the appendix B. 1 . 
Derivations of the beam filling factor Eq. (2), are described in the 
appendix B.2. 
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where D describes the diameter of the telescope, cught the speed 
of light, and ^ = 3600 ■ 180 'nT^ a conversion factor to get the 
telescope beam FWHM size in arcsec. For interferometric obser¬ 
vations, the user has to define the interferometric beam FWHM 
size directly. In contrast to single dish observations we assume 
a constant interferometric beam FWHM size for the whole fre¬ 
quency range. 

The term in Eq. (1) indicates the largest beam fill¬ 
ing factor of all core components of all molecules, that is, 
Cax = max{)7(6f"’")). 

m,c 

In general, the brightness temperature of radiation tempera¬ 
ture J{T, v) is defined as 


J(T,v)^ 


h V 


1 


kg g/iVAf _ 1 ■ 


(4) 


The expression 7cmb used in Eq. (1), describes the radiation 
temperature Eq. (4) of the cosmic background T^bg = 2.7 K, that 
is, 7cmb = J{Tcbg,v)- 


In Eq. (1), the expression represents the source func¬ 

tion and is according to Kirchhoff’s law of thermal radiation 
given by 


m,c / \ , m.c / \ 

e, (v) + e. (y) 

^ “ m,C/ \ , m.c / \ 

K, {V)+K^ (V) 

_ K, (V) /(Tex’ ,V)+K/ (V) J{T^’ , V) 


m.c / \ , 

K, iy) + K^ 

\ [r7-^(v)A7'er,v) + T7(y)y(7'7,v) 


d 


where ef^{v) and are the core and foreground coefficients 

for line and dust, respectively. Additionally, the optical depth 
is given by = f /(^’‘^(v)ds = a'"'’‘^(v) s. This assumes 

that molecules and dust are well mixed, that is, it would not 
be correct if the molecule exists only in part of the cloud, 
but the dust everywhere. In older versions, the background 
temperature could only be defined as the measured continuum 
offset, which corresponds to the beam-averaged continuum 
brightness temperature. At the same time, the dust, as agent of 
line attenuation, was described by column density and opacity. 
This is practical, because the observable Tbg is used, but does 
not constitute a self-consistent and fully physical description. 
Therefore, we now use optionally either a physical (y = 1) or 
phenomenological (y s 0) description of the background indi¬ 
cated by the Kronecker delta dyfl, that is, v) for 

j = 0. (Here, the phrase “background” means the “layer” with 
intensity which is located behind the core components, 

that is, the background of the core layers, see Eig. 1.) We note 
that, if y = 0, the definition of the dust temperature 7’J’‘^(v), 
Eq. (13), is superfluous. 


where the dust optical depth t“’^(v) takes the dust attenuation 
into account and is given by 


m,c 

(- 

' d,l-Qf 

Wref, 

Nm,c 

^Vref 


1 




4gas~dust 



(7) 


Here, describes the hydrogen column density, 

the dust mass opacity for a certain type of dust 
(Ossenkopf & Henning 1994) at the reference frequency 
Vref, and the spectral index of These parameters are 
defined by the user, see Sect. 2.1. In addition, Vref = 230 GHz 
indicates the reference frequency of the reference dust opacity 
'”^^2 describes the mass of a hydrogen molecule, and 

^gas-dust describes the dust to gas ratio and is set here to (1/100) 
(Hillebrand 1983). The equation is valid for dust and gas well 
mixed. 


The optical depth r”'’^(v) of all lines for each molecule m and 
component c is described as*^ 


“(v) = Xi 


87rv2 

X ^i'"’‘^’'(v) , 


-hv'/kB r; 


0 


( 8 ) 


where the sum with index t runs over all spectral line transitions 
of molecule m within the given frequency range. The Einstein 
Aui coefficient'^, the energy of the lower state £/, the upper 
state degeneracy g„, and the partition function Q {m, of 
molecule m are taken from the embedded SQLite3 database, 
described in Sect. 3. (Because the database usually does not 
describe the partition functions at the given excitation temper¬ 
ature , the value of Q {m, is computed from a linear 
interpolation. With the new catalog, extrapolation should not 
be necessary for most conditions encountered in molecular 
cores.) In addition, the values of the excitation temperatures 
and the column densities for the different components 
and molecules are taken from the user defined moffit file, see 
Sect. 2.1. 


In order to take broadening of lines caused by the thermal 
motion of the gas particles and micro-turbulence into account 
we assume in Eq. (8) a normalized Gaussian line profile, that is, 
(j)(v) dv - 1, for a spectral line t: 


0'"’^’'(v) 




cr" 


(v (v'^-8y^^sR)) 


(9) 


The source frequency for each component c of a molecule 
m is related to the user defined velocity offset taken 

from the aforementioned moffit file, by the following expression 


The total optical depth of each molecule m and com¬ 

ponent c is defined as the sum of the optical depths t”’^(v) of all 
lines of each molecule m and component c plus the dust optical 
depth T™’‘^(v), that is. 


total 


(v) = t;‘’"(v) + t7'(v), 


(6) 


_ '^oflFset , / 

CD - ' ^ 9 


^LSR 


Qight 


( 10 ) 


A derivation of the expression is given in the appendix B.3. 

The indices u and I represent upper and lower state of transition ?, 
respectively. 
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where v' indicates the frequency of transition f taken from the 
SQLiteS database mentioned above. Additionally, the standard 
deviation cr"'’‘' of the profile is defined by the velocity width 
described in the molht file for each component c of 
a molecule m\ 


Av”' 


rn.c.t _ 


( t , c m,cj\ 

^ +^^lsr) 


( 11 ) 


2 V2ln2 

The beam-averaged continuum background temperature 


/^“®(v) is parametrized as 


C'=(v) = 7'bg-( — 


slope 


+ -^CMB 


( 12 ) 


to allow the user to define the continuum contribution for each 
frequency range, individually. Here, Vmin indicates the lowest 
frequency of a given frequency range. Tbg and Tsiope, defined 
by the user, describe the background continuum temperature 
and the temperature slope, respectively. Here, the treatment 
of the dust is not entirely self-consistent. To amend that, we 
would need to define the sub-beam scale structure of the source, 
which we consider to be outside the scope of the current effort, 
although it is envisioned to provide this as an option in the future. 


In Eq. (1), the continuum contribution is described through 
the source function Eq. (5), by an effective dust tem¬ 
perature (through v)) for each component which is 

given by 

r;’^(v) = r:/(v) + 

( y \ d.slope 

— , (13) 

►^min / 

where and can be defined by the user for each 

component in the molht. If and T'”’’: are not dehned for 
a certain component, we assume = 7'™''^(v) for all com¬ 

ponents. Eor a physical (y s 1) description of the background 
intensity, see Eq. (5), the user can dehne the dust opacity, 
Eq. (7), and dust temperature, Eq. (13), for each component. 


Einally, the last term 7cmb in Eq. (1) describes the OEE 
position for single dish observations, where we have an intensity 
caused by the cosmic background 7cmb- For interferometric 
observations, the contribution of the cosmic background is 
hltered out and has to be subtracted as well. 


In contrast to core layers, foreground components may inter¬ 
act with each other, as shown in Eig. 2, where absorption takes 
places only, if the excitation temperature for the absorbing layer 
is lower than the temperature of the background. 

Hence, the solution of the radiative transfer equation for 
foreground layers can not be given in a form similar to Eq. (1). 
Eoreground components have to be considered in an iterative 
manner. The solution of the radiative transfer equation for fore¬ 
ground layers can be expressed as 

+ Cb'=(v) 

(v)«.,.=/ = 7 (0“’^=') (5'”'‘^=''(v) - r^r(v)«-,.=0-i)) 

X (l - (14) 



Fig. 2. Sketch of a distribution of core and foreground layers 
within the Gaussian shaped beam of the telescope (black ring). 
Here, we assume three different core components 2a, 2b, and 
2c located in a plane perpendicular to the line of sight which 
lies in front of the background layer 1 with intensity /^“‘^(v), 
see Eq. (12). The foreground layers 3 and 4 are located between 
the core layers and the telescope along the line of sight (black 
dashed line). Here, each component is described by different ex¬ 
citation temperatures, velocity offsets etc. indicated by different 
colors. The thickness of each layer is described indirectly by 
the total column density , see appendix B.3. Eor each fore¬ 
ground layer we assume a beam filling factor of one, i.e., each 
foreground layer covers the whole beam. 


where m indicates the index of the current molecule and i repre¬ 
sents an index running over all foreground components c of all 
molecules. Additionally, we assume that each foreground com¬ 
ponent covers the whole beam, that is, jj s 1 for all fore¬ 

ground layer. Thus, Eq. (14) simplifies to 

T^2"(y)>n.c=i = [5”'’-'(v)(i - + rx(v)e-"^o-'w] 

T^Z(y)>n,c=i = [5“’-'(y)(l O')) 

+ (15) 

where T^“‘^(v) describes the core spectrum, see Eq. (1), in¬ 
cluding the beam-averaged continuum background temperature 
/^™(v). Eor foreground lines the contribution by other compo¬ 
nents is considered by first calculating the contribution of core 
objects and then use this as new continuum for foreground lines 
reflecting the fact that cold foreground layers are often found 
in front of hotter emission sources. We assume, that the cosmic 
background describes together with the core components one 
end of a stack of layers. Additionally, the foreground com¬ 
ponents are located between this plane and the telescope, see 
Eig. 2. The total column density depends on the abundance 
of a certain molecule and on the thickness of a layer containing 
the molecule. The order of components along the line of sight is 
defined by the occurrence of a certain foreground component in 
the molht hie. 

By htting all species and their components at once, line 
blending and optical depth effects are taken into account. The 
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Rest Frequency (MHz) 

Fig. 3. Here, the myXCLASS function was used to calculate a 
spectrum for three different step sizes. 

modeling can be done simultaneously with isotopologues (and 
higher vibrational states) of a molecule assuming an isotopic 
ratio stored in the so-called iso ratio file, see Sect. 2.2. Here, 
all parameters are expected to be the same except the column 
density which is scaled by one over the isotopic ratio for each 
isotopologue. 

In order to correctly take instrumental resolution effects into 
account in comparing the modeled spectrum with observations 
myXCLASS integrates the calculated spectrum over each chan¬ 
nel. Thereby, myXCLASS assumes that the given frequencies v 
describe the center of each channel, respectively. This is partic¬ 
ularly important if the instrumental channel width is in the order 
of, or even larger, than intrinsic line widths. The resulting value 
is than given as 

1 r+T 

^mbCv) = — I T^b(v)dv, (16) 

AVc 

where Av^ represents the width of a channel. Due to the com¬ 
plexity of Eqs. (1) and (15) the integration in Eq. (16) can not 
be done analytically. Therefore, myXCLASS performs a piecewise 
integration of each component and channel using the trapezoidal 
rule and sums up the resulting values to get the final value used 
in Eq. (16). Internal resampling guarantees that even if the line 
width is smaller than the channel width the integration is done 
correctly, see Eig. 3. 

2.1. The moM file 

Within the molfit file the user defines both which molecules 
are taken into account and the number of components for 
each molecule. Additionally, the user has to define for each 
component the source size in arcsec (size), the excitation 
temperature Tex in K (T_ex), the column density Atot in cm“^ 
(N_tot), the velocity width (EWHM) Av in km s“' (V_width), 
the velocity offset in km s ' (V_o£f), which is via Eq. (10) 
connected to the source velocity Vlsr, and the flag (CFFlag) 
indicating if a component is considered for core c or foreground 
f. The definition of the source size for an foreground 


component will be ignored, because we assume that all fore¬ 
ground layers covers the whole beam. So, the definition of this 
parameter is not necessary. 

Example of a molfit file: 


% Number of molecules = 2 


% size: 

T_ex: 

N_tot: 

V_width: 

V_off: 

CFFlag 

CS;v=®; 

48.47 

3 

8®.®® 

3.91E+17 

2.86 

-20.56 

c 

21.8® 

51.®3 

6.96E+17 

8.®7 

30.68 

c 

81.7® 

68.11 

1.46E+17 

5.16 

-10.12 

c 

HCS-h;v=i 

®: 2 
15®.®® 

1. IOEh-18 

5.0® 

-0.15 

£ 


20®.®® 

2.20E+17 

3.1® 

-2.15 

f 


The definition of parameters for a molecule starts with a line 
describing the name of the molecule, which must be identical to 
the name of the molecule included in the database, see Sect. 3, 
followed by the number of components N for this molecule. 
The following N lines describe the parameters for each 
components, separately. Generally, all parameters have to be 
separated by blanks, comments are marked with the % character. 

In order to define the dust temperature for each component, 
the molfit file has to contain two additional columns between 
columns V_off and CFFlag, describing the parameters 
""d Eq- (13), respectively. 

Additionally, the myXCLASS program allows to define a hy¬ 
drogen column density Nh (in cm“^), dust mass opacity A:v„f (in 
cm^ g *), and the spectral index fi for each component or glob¬ 
ally, that is, —> Ah, a"™’, —> and —> (3. In order 

to define these parameters individually for each component, the 
molfit file has to contain additional columns on the left side of 
column CFFlag. 

Eor globally defined dust parameters, the myXCLASS program 
assumes that all core components do not contain dust except the 
core component with the largest beam filling factor, see Eq. (2). 
This avoids an overestimation of the dust contribution caused by 
the overlap of the core components. 

2.2. The iso ratio file 

As mentioned above, the myXCLASS program offers the possibil¬ 
ity to define isotopologues of a molecule and their abundance ra¬ 
tios to reduce the number of input parameters. Eor that purpose 
the user has to create an iso ratio file which consists of three 
columns, where the first two columns indicates the isotopologue 
and the corresponding molecule, respectively. The third column 
defines the assumed isotopic ratio. The columns are separated by 
blanks and comments are marked with the % character. 

Example of an iso ratio file: 


% isotopologue: molecule: ratio: 
S-33-02;v=®; S02;v=®: 3®.® 
S-34-02:v=®: S02;v=®: 5.® 
S02;v2=l; S02;v=®: 1.® 
S00-17:v=®; S02;v=®: 75®.® 
S00-18;v=®; S02;v=®: 5®®.® 


Here, we assume that the abundance of molecule 
"S02 ;v=®; " is 30 times higher than the abundance of its iso¬ 
topologue "S-33-02;v=0;". 
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Fig. 4. Here, the myXCLASS function was used to model HIFI 
data of Sgr B2m (black) using SO 2 (with three different compo¬ 
nents), SO (with one component), and HNO (with one compo¬ 
nent). The intensities of each component are shown in the bottom 
half. 


2.3. The myXCLASS function 

The myXCLASS function calculates a synthetic spectrum for a 
user defined frequency range, see Fig. 4. The function returns 
the calculated myXCLASS spectrum, a list of all transition fre¬ 
quencies within the defined range, and the intensities and opti¬ 
cal depths for each molecule and component as python arrays, 
which can be used for further analysis. An example call of the 
myXCLASS function is given in appendix A.2. 


3. Database 

The XCLASS interface contains a SQLiteS database in¬ 
cluding spectroscopic data from CDMS (Muller ef al. 2001; 
Muller et al. 2005) and JPL (Pickett et al. 1998). Data can be 
retrieved and updated via VAMDC, which is an interoperable 
e-Infrastructure for the exchange of molecular and atomic data 
between different databases. In principle, all databases which 
support the VAMDC standard and which contain the required 
spectroscopic information can be accessed via myXCLASS. A 
local database was chosen to be able to run the program without 
network access. CDMS also provides the latest version of the 
SQLiteS database, which is updated when new entries to the 
CDMS or JPL are made. In addition, the XCLASS interface 
allows in principle to add private entries (e.g., unpublished 
new spectroscopic data) to the database. New functions which 
facilitate adding private entries will be included in one of the 
next releases of the XCLASS package. 

The spectroscopic data is stored in two tables: The ta¬ 
ble partitionfunctions contains the partition functions 
2(m, Tex) for more than 1000 molecules for 110 different 
temperatures between 1.072 and 1000 K, which widely extends 
the standard range of the CDMS/JPL entries. This feature was 
deemed desirable because particularly for foreground lines with 
an excitation temperature Tex of 2.7 K, extrapolation from the 


canonical values (lowest calculated temperature of the partition 
function 9.75 K) could be wrong in either direction by large 
factors. Additonally, the table transitions contains more than 
5 415 000 transitions and includes the frequency v' (in MHz) 
with its uncertainty, the Einstein A„/ coefficient (in i“'), the 
upper state degeneracy g,,, the energy of the lower state £; (in 
K), and the quantum numbers for each transition. 

Entries in the tables are related by the species name, which are 
different from the names used in the Splatalogue catalog, for 
example, we use "HC-13-N; v=0; " instead of "HC(13)N v=0" 
(Splatalogue). In one of the next releases the Splatalogue names 
will be usable as well. 

Furthermore, the XCLASS interface provides functions to 
update and manage the myXCLASS database: The function 
UpdateDatabase updates the existing database file by down¬ 
loading the latest changes from the databases via VAMDC or 
by downloading a complete new database file from the CDMS 
server. The first option adds new entries to the tables and over¬ 
writes existing entries with new data. Private entries remain 
unchanged. The second option is much faster but removes all 
private entries by overwriting the database file with the most 
recent file from the CDMS server. In addition the function 
DatabaseQuery offers the possibility to send an user defined 
SQL query string to the database, for example, to add private 
entries or query about species names, content etc. Also, the in¬ 
terface contains functions ListDatabase, see appendix A.l, 
and GetTransitions, see Fig. 5, to display information about 
transitions within an user defined range. The GetTransitions 
function allows the user to select a region on the displayed spec¬ 
trum, and prints out all the lines from the catalog around the 
region selected. Additionally, the function offers the possibility 
to set filters in frequency and energy range. 

4. MAG IX 

So far, we have discussed how to use XCLASS to produce syn¬ 
thetic spectra. In fact, what the user usually desires is a de¬ 
scription of data that provides a fit to a set of observational 
data. In many packages, this is done by manually changing 
the parameters until the fit looks good. In XCLASS, we use 
an optimizer that provides the parameter set for a given model 
which best matches the data. Due to the large number of in¬ 
put parameters required by the myXCLASS program it is essen¬ 
tial to use a powerful optimization package to achieve a good 
description of observational data by the myXCLASS program. 
Therefore, the XCLASS interface contains the MAGIX package 
(Moller ef a/. 2013) which is a model optimizer providing an 
interface between existing codes and an iterating engine. The 
package attempts to minimize deviations of the model results 
from observational data, constraining the values of the model 
parameters and providing corresponding error estimates. 

MAGIX offers the possibility to model physical and chem¬ 
ical data using an arbitrary external model program, not only 
myXCLASS, explore their parameter space and find the set of pa¬ 
rameter values that best fits observational data. The goodness of 
a fit is described by the distribution which is a function of 
relative quadratic differences between observational and model 
values. MAGIX is able to explore the landscape of thex^ function 
without the knowledge of starting values. Additionally, MAGIX 
can calculate probabilities for the occurrence of minima in the 
X^ distribution and give information about confidence intervals 
for the parameters. The package provides optimization through 
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Fig. 5. Graphical user interface (GUI) of the GetTransitions 
function. Here, the user has to define a frequency range by se¬ 
lecting a single frequency in the shown spectrum. This frequency 
represents the central frequency of the range. The width of the 
range is given by an user defined parameter. XCLASS prints 
out informations (name of molecule, transition frequency, un¬ 
certainty of transition frequency, Einstein A coefficient, quantum 
number for lower and upper state, and energy of lower state) of 
all transitions in the selected molecule list located in the given 
frequency range, see screen output of ListDatabase described 
in appendix A.l. 


one of the following algorithms or via a combination of several 
of them (an algorithm chain or tree, see Fig. 6): the Levenberg- 
Marquardt (conjugate gradient) method, which is fast, but can 
get stuck in local minima, simulated annealing, which is slower, 
but more robust against local minima. Other global optimization 
algorithms, such as bees, genetic, particle swarm optimization, 
nested sampling, or interval nested sampling algorithms are in¬ 
cluded as well for exploring the solution landscape, checking 
for the existence of multiple solutions, and giving confidence 
ranges. Additionally, MGIX provides an interface to make sev¬ 
eral algorithms included in the scipy'® package available. 

Using an algorithm chain or tree offers the possibility to 
send the results of the optimization process performed by one 
algorithm to another optimization loop through some different 
algorithm. The simulated annealing as well as the Levenberg- 
Marquardt algorithm require starting values of the parameters 
that are not too far from the final solution, otherwise the algo¬ 
rithm might stray so far away from the solution that it never 
converges, that is, the user has to find an acceptable fit by hand 
before applying these algorithms produces useful results. Often, 
the location of the minimum can be guessed with sufficient ac¬ 
curacy to give good starting values, but sometimes one is com¬ 
pletely in the dark. Using an algorithm chain or tree, the user can 
first apply one of the swarm algorithms, for example, the bees 
or interval nested sampling algorithm, to determine the start¬ 
ing values for the subsequent local optimization algorithm using 
simulated annealing or the Levenberg-Marquardt algorithm, see 
Fig. 6a). But MAGIX does not only allow one to use the best but 
also the second best etc. result of a swarm algorithm as start¬ 
ing values for other algorithms (algorithm tree), see Fig. 6b). 
Therefore, MAGIX is able to find multiple minima of models. 


http://www.scipy.org/ 



Chain Tree 


Fig. 6. Example of an a) algorithm chain and b) algorithm tree. 
In both cases, we start with the bees swarm algorithm and use 
the best (a) and the three best results (b) as starting values for 
the Levenberg-Marquardt local optimization algorithm (LM), re¬ 
spectively. Afterwards, we apply the error estimation algorithm 
(Error) to the result(s) of the Levenberg-Marquardt algorithm. In 
case of an algorithm chain we get one final result (Rl), whereas 
we get three different results (Rl, R2, R3) for an algorithm tree. 


Since this can, depending on the data set and the number 
of free parameters, be quite computing intensive, OpenMP and 
OpenMPI parallelization for simultaneous evaluations of the 
external program is available for most of the algorithms. In 
addition to the normal MAGIX package which can be used 
with an arbitrary external model program, the XCLASS in¬ 
terface contains a further MAGIX version optimized for the 
usage with the myXCLASS program. This optimized MAGIX 
version is used by the myXCLASSFit and myXCLASSMapFit 
(myXCLASSMapRedoFit) functions and we are exploring the 
use of graphics processing units (GPUs). 

Generally, MAGIX is controlled by different XML files: The 
so-called instance file includes the names, initial values (and 
ranges) for all model parameters. Additionally, the instance file 
indicates the parameters which should be optimized by MAGIX. 
In addition the experimental XML file containing settings for the 
import of observational data, that is, path(s) and name(s) of the 
data file(s), format(s), range(s) etc. Furthermore, the algorithm 
control file, defining the algorithm or algorithm sequence which 
should be used for the optimization together with the settings for 
each algorithm. 

The new MAGIX function offers the possibility to use MAGIX 
within CASA with a previously registered external model pro¬ 
gram. For that purpose, the user has to provide the XML files, de¬ 
scribed above, and passes their path and file names to the MAGIX 
function, see A.3. 

In general, XCLASS creates job directories for some func¬ 
tions (MAGIX, myXCLASSFit, and myXCLASSMapFit) where all 
log and output files are stored. The job directories are created in 
the XCLASS run directory defined by the environment variable 
myXCLASSRunDirectory. The name of a job directory contains 
the date and time of the function execution followed by a unique 
job number. 
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5. Other functions 

The XCLASS interface for CASA includes two new functions 
(myXCLASSFit and myXCLASSMapFit) which provide a simpli¬ 
fied interface for MAGIX using the myXCLASS program. 

The myXCLASSFit function offers the possibility to fit mul¬ 
tiple ranges in multiple files from multiple telescopes simultane¬ 
ously using MGIX in conjunction with the myXCLASS program. 
It starts MAGIX using the Levenberg-Marquardt algorithm to op¬ 
timize the input parameters defined in an extended molfit file to 
achieve a good description of the observational data. The user 
has to specify the maximum number of iterations, the path and 
file name of the observational data file (or of a MAGIX exper¬ 
imental XML file, see Sect. 4), and the path and name of an 
extended molfit file. In contrast to the molfit file described in 
Sect. 2.1, the extended molfit file, see appendix A.4, contains 
one additional column on the left side of each parameter (except 
for the core/foreground flag) of each component. The additional 
column defines the range for each parameter. The definition of 
ranges guarantees that the fits does not run out of the defined 
parameter ranges. The limits, that is, the lowest and highest al¬ 
lowed values, for the parameters source size c, column density 
Atot, and hydrogen column density Nh are determined by 

lower limit = Iparameter value] / (range value) 

upper limit = Iparameter value] x (range value). 

The limits for the other parameters are calculated by simply 
adding or subtracting the value of the additional column to or 
from the value itself. All calculated lower limits which are less 
than zero are set to zero (except for velocity offset). If the range 
of a certain parameter is set to zero, the corresponding parame¬ 
ter is kept constant and is not optimized by the myXCLASSFit 
function. Additionally, the myXCLASSFit function offers the 
possibility to fit the iso ratios in the same fit process as well. For 
that purpose, the user has to add two additional columns in the 
iso ratio file on the right side of the third column. The fourth and 
fifth column defines the lower and upper limit for each iso ratio, 
respectively. If the lower limit is equal to the upper limit or if the 
lower limit is higher than the upper limit, the corresponding iso 
ratio is not optimized by the myXCLASSFit function. In general, 
the fit procedure stops, if the maximum number of iterations is 
reached, or if drops below 10 '^ (or an user defined value). 
Additionally, the myXCLASSFit function offers the possibility 
to use another algorithm or an algorithm chain or tree by 
defining the path and name of a MAGIX algorithm XML file, 
see Sect. 4. Furthermore, the function offers the possibility to 
fit multiple frequency ranges simultaneously by using a MAGIX 
experimental XML file, see Sect. 4. Finally, the myXCLASSFit 
function returns the optimized molfit file from the best fit and 
the corresponding modeled spectra as python arrays. (Here, the 
phrase “best fit” is connected to the application of an algorithm 
chain or tree: Using an algorithm chain or tree offers the 
possibility to use more than one algorithm to fit the data. Each 
application of an algorithm is connected with a fit represented 
by a certain value. The myXCLASSFit function reads in the 

values of all optimization loops of all applied algorithm and 
determines the lowest value. The result of the “best fit” is 
than given by the parameter vector which belongs to the lowest 
X^ value.) 

In addition to the myXCLASSFit function which is useful 
to fit single spectra, the XCLASS interface for CASA contains 
the myXCLASSMapFit function which fits one or more complete 
(FITS) data cubes. Here, the myXCLASSMapFit function reads 


Example of myXCLASSMapFit 
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Fig. 7. Example of parameter maps created by the 
myXCLASSMapFit function. Here, the K-ladder structure of the 
CH 3 CN( 12 - 11 ) transition (plus isotopologues) in G75.78-1-0.34 
was fitted at around 220 GHz (observed with the SMA). The 
temperature gradient is expected as also seen in different NH 3 
lines (Sanchez-Monge 2011; Sanchez-Monge et al. 2013). 


in the data cube(s), extracts the spectra for each pixel and fits 
these spectra separately using the Levenberg-Marquardt (or 
another) algorithm (chain or tree). The optimization procedure 
for each pixel stops, if the maximum number of iterations is 
reached, or if the x^ value drops below 10 “^ (or below an user 
defined value). In analogous to the myXCLASSFit function 
the fit parameters and their ranges for the myXCLASSMapFit 
function are defined in an extended molfit (and iso ratio) file 
as well. Additionally, the myXCLASSMapFit function offers 
the possibility to limit the fit to certain frequency ranges of 
a spectrum and to an user defined region of the cube(s). In 
order to reduce the computation effort, the user can exclude 
pixels by defining a threshold for the min. intensity of a pixel, 
that is, if the maximum intensity of a pixel is below the user 
defined threshold, the corresponding pixel is not fitted by the the 
myXCLASSMapFit function. In general, the myXCLASSMapFit 
function assumes, that the first three axes of the data cube(s) 
describe the right ascension, the declination, and the frequency, 
respectively. If the frequencies are not given in MHz, the EITS 
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header has to contain the definition of the CUNIT3 command 
which defines the unit of the frequency axis. If more than one 
data cube is specified in the XML file, the data cubes must 
describe the same map, that is, covered spatial area and the 
resolution have to be identical. The different data cubes are 
allowed to differ only in the frequency axis. At the end of the 
whole fit procedure, the myXCLASSMapFit function creates 
FITS images for each free parameter of the best fit, where each 
pixel corresponds to the value of the optimized parameter taken 
from the best fit for this pixel, see Fig. 7. Furthermore, the 
myXCLASSMapFit function creates FITS cubes for each fitted 
data cube, where each pixel contains the modeled spectrum. 
Finally, the myXCLASSMapFit function creates one FITS image 
which describes the quality of the fit for each pixel. Here, each 
pixel corresponds to the value of the best fit for this pixel. 
Applications of this are temperature maps, but also first and 
second moment maps, which are based on the simultaneous 
fitting of many lines, and are fairly robust against line confusion 
and blending of single lines. 

In order to improve the results of a previous application of 
the myXCLASSMapFit function the XCLASS interface contains 
the myXCLASSMapRedoFit function which offers the possibil¬ 
ity to redo one or more so-called pixel fits of a previous myX¬ 
CLASSMapFit run. The function performs fits for the selected 
pixels and recreates the different parameter maps using the new 
optimized parameter values. 


6. Example application 

The XCLASS interface for CASA was used to redo the 
analysis of the spectral line survey of the high-mass star¬ 
forming region Orion KL from 325 to 360 GHz as reported by 
(Schilke et al. 1997), which was done without using the myX- 
CLASS program nor optimization packages like MAGIX. 

In our fit, shown in Fig. 8 , we used the following 25 
molecules: C 2 H 3 CN, C 2 H 5 CN, CCH, CH 3 CCH, CH 3 CN, 
CH 3 OCH 3 , CH 3 OCHO, CH 3 OH, CN, CO, CS, H 2 CO, 
H 2 CS, H 2 O, HCCCN, HCN, HCO+, HCS+, HDO, HNCO, 
NO, OCS, SiO, SO 2 , and SO. The different molecules 
were described by 44 components in total, see Tab 1. 
Additionally, we have taken 241 isotopologues and higher 
vibrated states into account where we assumed the fol¬ 
lowing ratios: [^^C]/['^C] = 45 (Crockett ef al. 2014a), 

[''’0]/['*0] = 250 (Crockett effl/. 2014a), ['®0]/['’0] = 2625 
(TerceroefaL 2010), [^2S]/[33S] = 75 (Crockett ef aZ. 2014a), 
[32S]/[34S] = 22.5 (Schilke efflZ. 1997), ['4 n]/[' 5N] = 234 
(Crockett ef flZ. 2014a), pSi]/pSi] = 20 (Schilke ef aZ. 1997), 
and [^^Si]/[^°Si] = 30 (Anders & Grevesse 1989). In order to 
reduce the computational effort, we assumed a ratio of one for 
all vibrational excited transitions, that is, we used the same 
column densities as for the corresponding molecules. 

Additionally, we were able to (partly) identify 35 of the 57 
unidentified lines (U-lines) reported by (Schilke ef aZ. 1997), 
see Tab. 2. Here, we marked an unidentified line as (partly) 
identified, if the modeled integrated intensity of a component 
covers at least 10 % of the integrated line intensity 

(/cZv). We note that, the majority of the unidentified 
lines contain contributions from more than one component. 
Contributions of more than 100 % are caused by inaccuracies 
in the fit. For 19 U-lines the modeled spectrum describes less 
than 50 % of the integrated line intensity. Nevertheless, we are 



Rest Freauencv tGHzl 



337.50 337.82 338.15 338.47 338.80 



Rest Frequency (GFIz) 


Fig. 8. Result of the spectral line survey of the high-mass star¬ 
forming region Orion KL between 325 to 360 GHz using the 
myXCLASSFit function. Here, we used 24 molecules in conjunc¬ 
tion with 245 isotopologues. In panel a) the whole spectral line 
survey (back) is shown together with the result of myXCLASSFit 
function (red). Panel b), c) and d) show fit examples around 
332.5 GHz, 338.15 GHz, and 358.0 GHz, respectively. 


able to (partly) describe all unidentified lines, except the line 
at 348358 MHz, with an integrated line intensity of more than 
10 K km s *. The improved identification may caused by an 
updated and enhanced line catalog. 

In the following we will compare our fit results, given in 
Tab. 1, with those reported by (Schilke et al. 1997). Our calcu¬ 
lated errors describe the 2 cr credibility interval for each parame- 
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Table 1. Fit results of the molfit file parameters used for the Orion-KL model. Here, we describe only those molecules, which show 
at least three transitions within the covered frequency range. 


Molecule 

(arcsec) 

T5' (K) 

N]": (cm--) 

A v'"-'' (km s ') 

v;- (kms-‘) 

CHsCN 

0+37.6 

•^-2.0 

245+24.1 

^^*^-18.2 

3T7I4)™ 

10+4.0 

iU-22 

9+1.8 

'-1.3 


- 

[445:36] 

[1.7(15) 

- 

- 

CH3CCH 

70+9.8 

' -69.3 

91+26.3 

^^-14.5 

4.2 (16) :]::;]]> 

6+2.0 

^-6.1 

10+14 

^^-1.0 


- 

[65!^] 

[4.3(15) :];]<;|] 

- 

- 

HCCCN 

0 1+37.0 
° -13.6 

99+24.2 

'^-15.1 

i-o(15)!L1;:j 

11+3.5 

^^-3.1 

c+1.9 

*^-1.4 


9+35.8 

^-1.3 

330:f,° 

5.9(14):];];;® 

21+4.1 

^^-2.9 

4+1.7 

^-1.4 


- 

[225!™] 

[5.2(15):];];;®] 

- 

- 

CCH 

09 1 +29.0 
^ ' ^-326.1 

0 1 +265.5 

1-29.7 

1.6(11):-;*** 

0+4.3 

*^-3.0 

lO+ii.o 

^'^-10.3 



330:^5° 

7.9(14):];];;® 

0+2.4 

■^-2.7 

10+2.4 

^^-1.9 

OCS 


65+25.8 

^^^-13.7 

1.8(16):];];;]; 

^^-3.2 

9+1.6 

-1.3 


261 + 125.2 
^^^-260.0 


1-1(15):];];;® 

5+3.5 

*^-3.4 

0+2.2 

°-1.7 


- 

[83!“] 

[1-8(16) :;:];;6;] 

- 

- 

HNCO 

Q'T9+57.7 

-’-’^-331.0 

160+25.8 

1^11-20.7 

4 7 (14) 

'^^^2-1.3(14) 

6+3.6 
^-3 3 

9+2.0 

-1.3 


20+29.8 

^'^-16.7 

16+25.8 

^^-14.9 

7-4(16):;;];;]; 

^^-3.3 

9+1.8 

^-1.4 


- 

[i6o:]4] 

[1.3(15):];];;®] 

- 

- 


- 

[36:6] 

[5-7(15):];];;®] 

- 

- 


- 

[277:|«0] 

[1-3(15):®];;®] 

- 

- 


- 

[i9o:«] 

[1-2(15):®];;®] 

- 

- 

CN 

42+27.9 
^ -22.3 

40+21.9 

^*^-18.1 

1 7 (14)+^-2(‘5) 

A- / VA^;_3 703) 

4+4.0 

^-2.3 

10+1.7 
^ -2.1 


j 0+30.4 
^^-9.5 


1-3(15):];];® 

22^^-'^ 

^^-3.5 

0+2.7 

^-1.1 

H2CO 

249+2.2 

^^^-248.0 

i5o:3« 

8-4(14):];];;]; 

94+3.8 

^^-3.5 

0+1.9 

°-1.3 


97C+112.0 
-274.0 

130+23.1 

-24.6 

1.8(14):-;;® 

7+3.5 

^-3.5 

Q+2.0 

^-1.3 


- 

[Hl^CO: 230:l“] 

[h:3co:6.3(14):]];;:;] 

- 

- 

HDO 

^'^'^-75.4 

11+171.1 

^^-9.6 

5 0 (10) 

VAu; _5 ono) 

4+17.0 

^-4.4 

1 8+16 
^°-0.7 


47+12.4 
^ -45.5 

270:*’,3, 

5-7(15):®];;]; 

0+4.9 

^-2.7 

9+1.4 

-1.5 


- 

[272:“°] 

[0-9(16):;;];;®] 

- 

- 

H2CS 

ivj*j-]9 3 

105+26.7 

^'^‘^-14.7 

2.1(14):-;;® 

6+3.5 

^-3.4 

0+1.9 

°-1.3 


- 

[93:|°] 

[8-5(14):];];;:;] 

- 

- 

NO 

qo+299.8 

■^°-36.9 

47+49.2 

^'-15.0 

0 0 /1 91 +6.6(14) 

vAZ,; _3 302) 

10+1.2 

^'^-10.4 

5+3.2 

*^-0.2 


9Q+26.5 

^^-24.4 

196+23.9 

-21.5 

1 2 

A-^ vA / ; _3 

22^^-'^ 

^^-3.4 

0+1.3 

°-0.9 

SO 

jO+28.0 

^ -9.5 

^*^-14.7 

7 1(16)+9.3!l9i 
/•A VAU)_3mg+ 

^^^-2.1 

0+3.0 

^-0.9 


^+30.2 

^^-4.0 

264+20.5 

^^^-20.5 

2-9(16):];];;:; 

15+3.5 

^^^-2.4 

0+1-8 

^-1.9 


- 

[27:]] 

[6-3(16):®];;®] 

- 

- 

CH3OCHO 

0+38.8 

•^-2.3 

1 1 9+24.0 
^ ^ -20.6 

Q 0 /I +6.6(18) 

VA^; -3.4(15) 

c+3.9 

*^-2.4 

Q + 1.0 
^-1.1 

CH3OCH3 

O+40.0 

•^-2.0 

[98:3] 

136+24.6 

IJU 20.4 

[1-5(16):®;;;®] 

1 0 (16) 

A.AJ VAU)_J 

0+3.9 

■^-2.5 

q+1.4 

^-0.9 


94+25.9 

^^-23.2 

94+22.5 

^^-22.3 

2.8(15)::°®;:® 

6+3.9 

^-2.5 

0+1.3 

^-1.0 


- 

[89:6] 

[1-8(16):®];;®] 

- 

- 

C2H5CN 

90+15.0 

^^-35.5 

127+23.3 

^^^-21.6 

9-9(14):]:;;;:; 

^"^-3.1 

5+1.7 

*^-1.4 


- 

[99:3] 

[1-3(16):®];;®] 

- 

- 

C2H3CN 

15+29.1 

^*^-13.6 

174+31.4 
^ ^^-16.1 

8.7(13):-;;® 

4+1.0 

^-0.7 

6+1.0 

^-0.9 


1 0+29.2 
^^-17.2 

^^^-21.5 

3 5 (14) 

VA-t; _7 6031 

12+i-O 

^^-0.7 

0+0.9 

“^-1.2 

CH 3 OH 

1 9+27.6 
^^-10.5 

[96:6] 

9^+25.8 

^*^-14.4 

[8.2(14):;;];;:]] 

9 2 (15) 

VA^; -2.8(15) 

5+3.7 

*^-2.6 

q+1.8 

^-1.4 


0+37.1 

*^-2.1 

189+26.0 

iO“-20.0 

4 5 (16) 

‘A.J VAU) _2.1(16) 

0+3.8 

°-2.5 

9+1.5 

^-1.4 


4+36.9 

^-2.5 

161+28.1 

^^^-18.3 

1-3(16)::;];;® 

9+3.9 

^-2.4 

0+2.1 

^-1.8 


- 

[188:]] 

[7-0(16):®];;®] 

- 

- 

SO 2 

O+40.6 

*^-1.9 

140+29.7 

^^^-16.5 

2 3 (16) 

VAu; _9 3(15) 

6+3.9 

^-2.4 

c+1.0 

*^-1.1 


q+29.3 

^-7.8 

1 90+29.9 
^ ^'^-18.9 

^•A -2.1(16) 

0 1 +4.6 
■^-^-2.8 

(2+0.8 

^^-1.1 


16+29.4 

^^-14.6 

90+25.4 

^"^-14.9 

5.4(16):];];“; 

9 1 +4.9 
^^-2.5 

(2+0.8 

^^-1.1 


- 

[124:]] 

[7-7(16):®];;®] 

- 

- 


Notes. The column densities are given on a main beam brightness line temperature scale, i.e., each column density (exponents are described by 
round brackets, e.g., 1.4(16) = 1.4 x lO'*) is multiplied with the beam filling factor of the corresponding component, to make them comparable 
to the results given by (Schilke et al. 1997), which were indicated here within squared brackets. Additionally, the errors of the fit parameters 
(indicated by the sub-(left) and superscript (right) values) were determined using the error estimation algorithm included in the MAGIX package. 
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ter, respectively. Large asymmetric error ranges for some column 
densities (e.g., SO) are caused by the fact that some lines start 
to become optically thick which prevents an accurate parameter 
estimation. For these lines we can only determine lower limits 
of the corresponding column densities. 

In contrast to (Schilke et al. 1997) we were not able to iden¬ 
tify NH 2 CHO, C 2 H 5 OH, and HCOOH because of their weak 
lines which do not allow a distinct estimation of their contribu¬ 
tion to the given spectrum. But we took CCH, CN, and NO into 
account. 

Our derived excitation temperatures and beam averaged col¬ 
umn densities for H 2 CS and CH 3 OCHO agree (within the given 
error ranges) with Schilke’s result. In order to achieve a good 
description for HCCCN, HDO, and OCS we used two compo¬ 
nents, respectively. For HCCCN and HDO we find that the sec¬ 
ond components show similar temperatures and column densi¬ 
ties compared to Schilke’s results. Due to the large error ranges 
for OCS we are not able to identify the corresponding compo¬ 
nent for Schilke’s result although the temperatures as well as 
column densities are comparable within the given error ranges. 

Additionally, we used H 2 CO with two components instead 
of H^^CO with one component. Using the aforementioned iso- 
topologue ratio of ['^C]/['^C] = 45 we find that the temperature 
and column density for the first component of H 2 CO agree with 
those described by (Schilke et al. 1997). 

As shown in Tab. 1, we used three components to de¬ 
scribe the contribution of CH 3 OH and SO 2 . Schilke’s result 
for CH 3 OH corresponds to our second component whereas 
the excitation temperature and beam averaged column density 
for the first component of SO 2 agree with those reported by 
(Schilke et al. 1997). Furthermore, the described temperatures 
and column densities (after scaling with the isotopologue ratios, 
mentioned above) for ^‘*S 02 match our results nicely. For ^^S 02 
the derived excitation temperature do not coincide with Schilke’s 
result. 

In order to model the contribution of HNCO we used only 
two instead of four components. We see that the temperatures 
and column densities of Schilke’s first, third, and fourth compo¬ 
nent match our results for the first component. But, the column 
density of Schilke’s second component clearly does not agree 
with our column density for the second component although the 
excitation temperatures are comparable. 

Although we find similar results for the majority of 
molecules we derived different excitation temperatures for 
CH 3 CN, CH 3 CCH, SO, CH 3 OCH 3 , C 2 H 3 CN, and C 2 H 5 CN. In 
contrast to CH 3 CCH, SO, C 2 H 3 CN, C 2 H 5 CN, and CH 3 OCH 3 
where we find discrepancies with less than 24 K we find a strik¬ 
ing difference for CH 3 CN. Already (Schilke et al. 1999) noticed 
that the value of 445 K based on rotation diagram analysis was 
wrong due to a neglect of the line opacity. They presented a 
corrected fit using the opacities derived from CH 3 CN/CH 3 ^CN, 
yielding 160 K, which differs from the present value because 
(i) a [^^C]/['^C] ratio of 60 (instead of 45) was used then 
- which increases the correction due to opacity, and (ii) the 
correction was done by hand, not taking the individual lines as 
well into account as the present study. Not taking into account 
the opacities also underestimates the column density, by an 
order of magnitude in the present case. Additionally, we derived 
a column density for CH 3 CCH which is nearly an order of 
magnitude higher than Schilke’s result. 


Table 2. New identified lines. Here, an unidentified line is 
marked as identified, if the modeled integrated intensity of 
at least one molecule covers at least 10 %. The numbers in 
squared brackets describe the components of the corresponding 
molecules, respectively. 


V 

(MHz) 

j'TSchilke^y 

(K km s“*) 

(K km s“*) 

Molecule 

330715 

1.7 

3.3 (195.7 %) 

CHiOCHO,Vi 8 = 1 , [ 1 ] 

332789 

3.8 

1.0 (25.9%) 

C 2 H 3 CN, [2] 

333865 

8.9 

1.4(15.5 %) 

^^SO, [ 1 ] 

334140 

1.6 

1.5 (91.9%) 

CH30CH0,Vi8 = 1 , [ 1 ] 

335335 

6.2 

0.8 (12.7 %) 

C 2 H 5 CN, [1] 

335703 

4.6 

1.0 ( 22 . 2 %) 

CH 3 OH, [2] 

335840 

2.9 

0.8 (27.3 %) 

CH 3 OCHO, [1] 

337744 

4.2 

1.3 (31.6%) 

CH30H,Vi2=1, [ 2 ] 

337839 

3.3 

4.9 (149.7 %) 

HCCCN, V7 = l, [2] 



3.1 (94.2%) 

CH 3 OH, [2] 



0.8 (23.0 %) 

CH 3 OH, [3] 

339138 

3.3 

0.5 (14.0 %) 

‘ 3 CH 3 CN-A [1] 

339527 

8.6 

1.2(13.4%) 

CN, [2] 

340496 

15.9 

6.9 (43.5 %) 

C 2 H 5 CN, [1] 

340872 

4.6 

0.7 (14.4 %) 

^^SO, [ 1 ] 

341472 

10.9 

10.3 (94.9 %) 

C 2 H 5 CN, [1] 

341482 

11.7 

2.1 (18.1 %) 

C 2 H 5 CN, [1] 

341499 

0.4 

0.2 (54.1 %) 

CH30CH0,vi8 = 1 , [ 1 ] 

342129 

4.1 

0.6 (14.7 %) 

C 2 H 3 CN, [2] 



0.7 (16.2 %) 

CH30CH0,vi8=1, [ 1 ] 



0.4 (10.3 %) 

C 2 H 3 CN, [1] 

342290 

3.9 

3.7 (96.1 %) 

CH30CH0,vi8 = 1 , [ 1 ] 

342486 

4.7 

0.9 (20.0 %) 

CH30CH0,Vi8=1, [ 1 ] 

343202 

21.8 

13.1 (60.2%) 

C 2 H 5 CN, [1] 

343665 

7.8 

1.7 (21.4%) 

CH30CH0,Vi8 = 1 , [ 1 ] 

344773 


0.2 

[ 2 ] 

344788 


1.2 

[ 2 ] 



1.0 

CH 3 OCHO, [1] 



0.2 

^'^SOz, [3] 

344796 


2.5 

^'‘SOz, [ 2 ] 



0.8 

^'^SOz, [3] 

347446 

4.7 

2.2 (45.9 %) 

C 2 H 3 CN, [2] 



0.8 (17.8 %) 

CH30H,Vi2=1, [ 2 ] 



0.7 (15.0 %) 

C 2 H 5 CN, [1] 

348084 

2.4 

0.5 (22.7 %) 

^'‘SOz, [ 2 ] 



2.3 (95.5 %) 

CH30CH0,vi8 = 1 , [ 1 ] 

348373 


31.0 

SO 2 , [ 2 ] 



0.2 

SO‘^0, [3] 



0.2 

HCCCN,V 7 =3, [2] 



0.2 

SO'^^O, [ 2 ] 



0.4 

SO'^^O, [ 1 ] 



0.2 

CCH,V 2 = 2 , [ 2 ] 



0.3 

C 2 H 5 CN, [1] 

350170 

2.8 

2.3 (81.3 %) 

CH3CN,V8 = 1 , [ 1 ] 



1.5 (53.3 %) 

C 2 H 5 CN, [1] 

351490 

7.2 

1.1 (15.8%) 

“SOz, [3] 

351540 

10.8 

5.0 (45.9 %) 

HNCO, [1] 



3.1 (28.6%) 

“SOz, [3] 



1.3(11.9%) 

C 2 H 5 CN, [1] 

353166 

5.0 

1.2 (23.3 %) 

C 2 H 3 CN, [2] 

354129 

7.8 

5.6 (72.3 %) 

CH30H,Vi2=1, [ 2 ] 



1.2 (14.7 %) 

CH30H,vi2=l, [3] 

355851 

3.2 

0.4 (12.8 %) 

OC^^S, [ 1 ] 

356644 

5.3 

1.1 (21.4%) 

HCCCNv6=1,V7=1,[2] 

358356 

14.1 

1 . 6 ( 11.1 %) 

^'‘SOz, [ 2 ] 



1.9 (13.5 %) 

^'^SOz, [ 1 ] 



3.3 (23.5 %) 

C 2 H 5 CN, [1] 


For optically thin, unblended lines (like for H 2 S Notes. Due to inaccuracies in the fit we overestimated the integrated 

and CH 3 OCHO), the rotation diagram method gives hne intensities of some weak lines, 

the same results as the XCLASS fit, but requires more 
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effort through manual fitting of all lines. High opac¬ 
ities, while they can be corrected in rotation diagrams 
(Schilke ef fl/. 1997; Goldsmith & Langer 1999), this is very 
cumbersome to do by hand, and is much more robustly achieved 
through the use of XCLASS - in this data set, CH 3 CN and 
CH 3 OH are prominent examples. Lastly, the identification of 
species with weak features is only possible reliably in line-rich 
spectra if all the other species are modeled as well, as can be 
seen by the XCLASS fits reported in (Belloche ef aZ. 2008; 
Belloche ef aZ. 2009; Belloche ef aZ. 2013). All this together 
with the capability of producing credibility intervals demon¬ 
strates that the methods of XCLASS give much more accurate 
and reliable results than those achieved with rotational diagram 
methods. 


7. Conclusions 

We presented the XCLASS interface for CASA which provides 
powerful new functions for the CASA distribution for analyzing 
spectral surveys. The toolbox includes the myXCLASS program 
which models observational data by solving the radiative trans¬ 
fer equation for an isothermal object in one dimension (detec¬ 
tion equation). The myXCLASS program is able to model a huge 
number of molecules simultaneously, whereat the contribution 
of each molecule can be described by an arbitrary number of 
components. These can usually be distinguished by different ra¬ 
dial velocities and do not interact with each other radiatively but 
superimpose in the model. The myXCLASS program depends on 
a number of parameters, which are partially taken from an em¬ 
bedded SQLiteS database containing entries from CDMS and 
JPL through the VAMDC portal. The other input parameters are 
read from an user defined ASCII file containing parameters for 
each molecule and component. In order to achieve a reason¬ 
able description of observational data one has to partially op¬ 
timize the user defined parameters. For that purpose the toolbox 
contains an interface for the model optimizer package MAGIX, 
which helps to find the best description of the data using a cer¬ 
tain model, that is, finding the parameter set that most closely 
reproduces the data. Therefore, the toolbox provides two func¬ 
tions to model one or more single spectra or data cubes simulta¬ 
neously using the myXCLASS program in conjunction with one of 
the optimization algorithms (or combinations of them in an algo¬ 
rithm chain or tree) included in the MAGIX package. The combi¬ 
nation of myXCLASS and MAGIX prepare the way for an automatic 
line identification routine which is essential for a spectral sur¬ 
vey on large spectral cubes produced by ALMA. Additionally, 
the MAGIX package offers the possibility to apply more elaborate 
programs (such as LIME) to model astronomical data as well. 

Enhancements are being worked on a 3-d version for defin¬ 
ing sub-beam structure. Additionally, we work on an automatic 
line identification function for XCLASS, which will be de¬ 
scribed in a subsequent paper. 
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Appendix A: Examples 

A.1. ListDatabase function 
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This function reads in entries from the table transitions, see Sect. 3, located in the SQLiteS database file and prints out the 
contents (name of molecule, quantum number for lower and upper state, transition frequency, uncertainty of transition frequency, 
Einstein A coefficient, and energy of lower state) to the screen or file (defined by the input parameter OutputDevice). The user can 
limit the output by defining a minimum and maximum for the frequency (or for the lower energy) for the transitions or by selecting 
one or more molecule(s) by using the input parameter SelectMolecule. 

Input within CASA: 

FreqMin = 2Q88Q.Q 
FreqMax = 998888.8 
ElowMin = 188.8 
ElowMax = 2888.8 
SelectMolecule = ["C0:v=8:"] 

OutputDevice = "" 

Contents = ListDatabase() 


Screen output (shortened): 


Name: 

upper state: 

lower state: 

Ereq. [MHz]: 

Error Freq. [MHz]: 

Einstein A[s-1]: 

E_low [K]: 

C0:v=8; 

X;v=8:3=7; 

X;v=8:3=6; 

886651.88688 

5.8888e-83 

3.422e-85 

115.354 

C0:v=8; 

X:v=8:3=8; 

X:v=8:3=7; 

921799.78888 

5.8888e-83 

5.134e-85 

153.798 


A. 2. my XCLASS function 

The function calculates a spectrum based on a user defined molfit file described in Sect. 2.1. The function returns the modeled 
spectrum as well as the intensities and optical depths for each component. 

In the example described below, we use the molfit file CH30H__pure .molfit located in "demo/myXCLASS/" to create a spec¬ 
trum between 580102.0 and 580546.5 MHz with a step size of 0.5 MHz. Additionally, we use global definitions for Tbg, Tsiope, Nh, 
and /3 (indicated by t_back_£lag = T and nH_£lag = T). Furthermore, we use the iso ratio file isojnames. txt and set the 
diameter of the telescope to 3.5 m. In addition, we set the rest frequency to 8.8 MHz with a reference velocity of 0 km s '. The 
output parameter modeldata contains the modeled spectrum described by a python array with two columns, where the first column 
indicates the frequencies (in MHz) and the second column the corresponding intensities (in K). For a rest frequency unequal zero, 
the my XCLASS function adds a further column to parameter modeldata describing the velocities (in km s'). 

The output parameter log contains some informations about the calculation and parameter trans describes the transition 
frequencies of all molecules used in the molfit file within the given frequency range. The intensities and optical depths for each 
component and molecule is stored in parameter IntOptical. Finally, parameter jobDir indicates the path and name of the 
corresponding job directory containing all files produced during the calculation process. We note that the user is free to define other 
names for the output parameters. 

Input within CASA: 

FreqMin = 588182.8 
FreqMax = 588546.5 
FreqStep = 5.8888888888E-81 
TelescopeSize =3.5 
Inter_Flag = F 
t_back_£lag = T 
tBack = 1.1 

tslope = 8.8888888888E+88 
nH_£lag = T 
N_H = 3.8888888888E+24 
beta_dust =2.8 
kappa_1388 = 8.82 

MolfitsFileName = "demo/myXCLASS/CH30H_pure.molfit" 

iso_£lag = T 

IsoTableFileName = "demo/myXCLASS/iso_names.txt" 

RestFreq =8.8 
vLSR =8.8 

modeldata, log, TransEnergies, IntOptical, jobDir = myXCLASSO 
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A.3. MAGIX function 

The function starts the MAGIX program. The user has to define the paths and names of the registration, instance, experimental, and 
algorithm xml file, respectively. The parameter MAGIXOption indicates if informations are printed to screen (None) or not (quiet). 
The results, that is, the optimized instance file(s) etc., are stored in a job directory described by parameter JobDir, mentioned above. 

Input within CASA: 

MAGIXExpXML = "demo/MAGIX/TwoOscillators_Re£Fit_R.xml" 

MAGIXInstanceXML = "demo/MAGIX/parameters.xml" 

MAGIXFitXML = "demo/MAGIX/Levenberg-iy[arquardt_Parameters. xml" 

MAGIXRegXML = "Fit-Functions/Drude-Lorentz_conv/xml/" 

MAGIXRegXML += "Conventional_Drude-Lorentz.xml" 

MAGIXOption = "" 

JobDir = MAGIX() 


A.4. myXCLASSFit function 

This function provides a simplified interface for MAGIX using the myXCLASS program. The function returns the contents of the 
molfit file containing the parameters for the best fit and the corresponding model function values for each data point. 

Here we fit observational data stored in file bandlb.dat using the molfit file CH30H__old.molfit, see below, without 
isotopologues (indicated by iso_flag = F) and the Levenberg-Marquardt algorithm with maximum ten iterations. 

Input within CASA: 

Numberlteration = 18 

MolfitsFileName = "demo/myXCLASSFit/CH30H_old.molfit" 

experimentalData = "demo/myXCLASSFit/bandlb.dat" 

TelescopeSize =3.5 
Inter_Flag = F 
t_back_flag = T 
tBack = 1.1 

tslope = 8.8808888888F+88 
nH_flag = T 
N_H = 3.08888Q8888F+24 
beta_dust =2.8 
kappa_1388 =8.82 
iso_£lag = F 
IsoTableFileName = "" 

RestFreq =0.8 
vLSR =0.0 

newmolfit, modeldata, JobDir = myXCLASSFit() 


Extended molfit file CH30H__old.molfit used in the example described above: 


CH30H;v=0; 1 

% limit: s_size: 
10.00 1.153 


limit: T_rot: 

50.00 115.1219 


limit: N_tot: 

50.00 1.648e+17 


limit: V_width: 
10.00 7.463 


limit: V_off: CFFlag: 

0.00 -7.900 c 


A. 5. myXCLASSMapFit function 

The myXCLASSMapFit function uses the myXCLASS program in conjunction with MAGIX to fit complete data cubes instead of single 
spectra. In the example described below, we fit a part of the FITS cube Orion. methanol. cbc. contsub. image. fits described 
by the region file region__box__ds9.reg. (In order to fit the whole data cube, the parameter regionFileName has to be cleared, 
that is, regionFileName = "".) Here, we set the threshold to two, that is, we fit all pixel which have a maximum intensity of 
0.001 K (or higher). Additionally, we use the molfit file CH30H.molfit together with the iso ratio file iso_names.txt and set 
the diameter of the telescope to 3.5 m. Furthermore, we use the Fevenberg-Marquardt algorithm with maximum 20 iterations for 
each pixel fit. In order to use another algorithm (or algorithm combination) the parameter AlgorithmXMLFile has to describe the 
path and name of a MAGIX algorithm XMF file. The myXCLASSMapFit function returns the path and the name of the job directory 
described by parameter JobDir containing subdirectories for each pixel fit and the parameter FITS images. 


14 



T. Moller et at: XCLASS 


Input within CASA: 

Numberlteration =15 

AlgorithmXMLFile = "demo/myXCLASSMapFit/algorithm-settings.xml" 

MolfitsFileName = "demo/myXCLASSMapFit/CH30H_new.molfit" 

experimentalData = "demo/myXCLASSMapFit/Orion.methanol.cbc.contsub.image.fits" 
regionFileName = "demo/myXCLASSMapFit/ds9_phys.reg" 

FreqMin =0.8 
FreqMax =0.0 
TelescopeSize =3.5 
Inter_Flag = F 
t_back_flag = T 
tBack = 9.5000000000E-01 
tslope = O.OOOOOOOOOOE+00 
nH_£lag = T 
N_H = 3.OOOO0OOOOOE+24 
beta_dust = 1.4 
kappa_1300 = 0.02 
iso_£lag = E 
IsoTableFileName = "" 

RestFreq =0.0 
vLSR =0.0 
Threshold = 0.001 
UsePreviousResults = T 

clusterdef = "demo/myXCLASSMapFit/clusterdef.txt" 

JobDir = myXCLASSMapFitO 


A.6. myXCLASSMapRedoFit function 

This function offers the possibility to improve the result of a previous performed myXCLASSMapFit run using another molfit file 
and/or other algorithms. 

In the example given below, the myXCLASSMapRedoEit function will re-fit the pixels (83.201,-15.1345) and 
(83.111,-15.221) of the previous myXCLASSMapFit run with job number 1234 using the molfit file "CH30H.molfit". 

Input within CASA: 

JobNumber = 1234 

PixelList = [[83.201, -15.1345], [83.111, -15.221]] 

Numberlteration = 10 
AlgorithmXMLFile = "" 

MolfitsFileName = "demo/myXCLASSMapFit/CH30H.molfit" 
experimentalData = "" 

Threshold =0.0 

FreqMin =0.0 

FreqMax =0.0 

t_back_flag = T 

tBack = 9.5000000000E-01 

tslope = O.OOOOOOOOOOE+00 

nH_£lag = T 

N_H = 3.OOOO0OOOOOE+24 

beta_dust =2.0 

kappa_1300 = 0.02 

iso_£lag = T 

IsoTableFileName = "demo/myXCLASSMapFit/iso_names.txt" 

RestFreq =0.0 
vLSR =0.0 
myXCLASSMapRedoFit0 
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Appendix B: Derivations 

B.1. Detection equation 

In absence of scattering, the radiative transfer equation 

^ = -Ky{s)Iy + ey(s) (B.l) 

di 

describes the propagation of radiation which passes through a medium. During the propagation photons are absorbed and emitted 
indicated by the emission and absorption coefficients Cy and Ky, respectively. 

The optical depth Ty(s) which measures the distance in units of the mean free path, is given by 

ps'=s 

Ty(s) - I Kyds'. (B.2) 

By using the source function Sy - the radiative transfer equation can be expressed as 

^^Iy(T)+Sy(T). (B.3) 

dr 

Within the local thermodynamic equilibrium (LTE) the source function is described by Planck’s law for black body radiation By{T). 


Integrating Eq. (B.3) over t leads to 

Iy{s)^Iy(0)e-^''+ f Sy{T'y^-^^dT'^. (B.4) 

Jt'=0 

Assuming a constant source function, that is, constant emission and absorption coefficients through the medium, the transfer equa¬ 
tion can be written as 

I^{s) = ly(0) +Sy(l- ■ (B.5) 

Additionally, we have to take into account, that the different components may not cover the whole beam, that is, that the background 
behind a certain component might contribute as well. So, we have to extend Eq. (B.5) by introducing the beam filling factor 77 , 
Eq. (2), which describes the fraction of the beam covert by a component 

/y(s) = 77 [/v(0)e-"' + (1 - ey + (1 - 77 ) 7,(0). (B. 6 ) 

Here, the term r//y(0)e~^’' indicates the attenuated radiation from the background /v(0). The second term 775v (1 - e”’^') describes 
the self attenuated radiation emitted by a certain component. Einally, the last term (1 - 77 ) /v(0) represents the contribution of the 
background which is not covert by a component. 

In real observations, we do not measure absolute intensities but only differences of intensities, that is, we have to subtract the 
OEF-position for single dish observations from Eq. (B. 6 ) as well, where we have an intensity caused by the cosmic background 
■fcMB- So, we achieve 


/y(s) = 77 [/v(0) e +Sy(l-e ’■■')] -H (1 - 77 ) /y(0) - dcMB 

- 77 [(5, - 7,(0)) (1 - ey] + 7,(0) - 7 cmb. 


(B.7) 


B.2. Beam fiiiing factor 

The derivation of the beam filling factor (2) starts with the normalized two-dimensional Gaussian function 




1 


^27T [o-l + o-fj 




(B. 8 ) 


where cr^ and cr^ describe the variances and /r , and fiy the center along the x and y axis, respectively. 

Observing a Gaussian shaped extended source with a telescope is described by a convolution of two two-dimensional Gaussian 
functions; 

X OO ^00 

I gi(x- u,y - v,cr^^u Cy,! , gliu, V, a-^^2, o-yg,IJ.X.2,Byg) du dv 

00 \ J —00 




(7J.r,l+7'x.2--^)^ (7'y.l+7'y.2-^-v)^ 


(B.9) 
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Fig. B.l. a) A single two-dimensional Gaussian function, Eq. (B.8), with cTx - (Ty - 1.9 and — Av - 0 representing a Gaussian 
beam of a telescope, b) Cut through the Gaussian function described in a) at half height, c) A single two-dimensional Gaussian 
function, with cr^ = cTy - 1.1 and - 0.5 and fiy - 0.9 representing a Gaussian beam of a point source, d) Cut through the 
Gaussian function described in c) at half height, e) Convolved Gaussian of telescope and point source as given by Eq. (B.IO). f) Cut 
through the convolved Gaussian function described in e) at half height. 

Assuming that gi describes the telescope with pix,\ - f^y.i = 0 and that telescope and extended source are described by non-elliptical 
Gaussians, that is, (Tx,i = cr^.i = cri and o-x ,2 = o-y ,2 = o' 2 , Eq. (B.9) can be simplified to*^ 

e . (B.IO) 


(^1 * 82) 


27t(^ 


cr, + crt 
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Here indicates the convolution of two functions gi and g 2 . 
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The FWHM of the resulting Gaussian is given by 


FWHM = 2 V2 log 2 ^cr\ + cr: 


(B.ll) 


which describe an area of 


= TT ■ 2 log 2 ■ (crj + cri) 




(B.12) 


In the last line we used the relation between the variances (Ti ,2 and the user defined FWHM of telescope {0\ = 0,) and source size 
(02 = 0'”’^^) which is given by 

0; = 2 V2 log2 ■ cr,-. (B.13) 

The beam filling factor Eq. (2) is defined as ratios of areas 


%l,g 2 ” 


“ 4 “ 


iFWHM 


|.(02 + 02 ) (02 + 02 )’ 


(B.14) 


which is completely independent of the position and /r,, of the extended source within the telescope beam. 


If more extended sources are observed with the telescope, we have to convolve the already convolved Gaussian Eq. (B.IO) with 
a further two-dimensional Gaussian function with variance cr+3 = cr^, 3 s 0-3 and center and 3. We get 


((.gl * g2) * g3) = (gl * g2 * gs) = 


1 


cr\ + cr\+ afj 


2ti 

V 2 log 2 -y^ 


2 (., 2 +£r 2 +cr 2 ) 


(j\ + cr\ + CTj 


The FWHM of the resulting Gaussian is given by 

FWHM = 2 

which describe an area of 

= TT ■ 2 log2 ■ ( 0-2 -H 0-2 -H crl) 


(B.15) 


(B.16) 


. ^.(02 + 02 + 02 )_ 


(B.17) 


In the last line we used again the relation between the variances and FWHM, Eq. (B.13). Now, we can again define a beam filling 
factor similar to Eq. (2) and get 

^ nl ,n2 

„ _ _ 4 _ 4 _ _ 2 3 /Tj 1 o\ 

^gl,g2,g3 (B.lb) 


|. (02 + 02 + 02 ) (02 +02 + 02 )- 


When we observe more than three extended sources with a telescope we have to iteratively convolve the resulting Gaussian with a 
further two-dimensional Gaussian function and we can generalize Eq. (B.18) to 


^general — 




Zii>i of 


{lliOf) 0f + Zi>i0f’ 

where 0i = 6, describes the FWHM of the Gaussian shaped beam of the telescope, defined by the diffraction limit, Eq. (3). 


(B.19) 


B.3. Optical depth 

In order to derive Eq. (8) we consider a system which involves radiative transitions between a lower I and an upper u level only. As 
shown in Eig. B.2, the lower level has an energy £; and the upper level an energy £„ > With 

hv,j^E,-E, (B.20) 

describing the energy difference between these two levels we can express the emissivity due to spontaneous radiative decay as 

h V 

ei,u,v ^ ^nuAu,i(l)i,u{v), (B.21) 

4n 

where A„j describes the Einstein A-coefficient, or radiative decay rate for the transition from the lower I to the upper u level. The 
expression J gives the averaged time, that a quantum mechanical system can stay in level u before radiatively decaying to level /, 
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Fig.B.2. Transitions between the lower I and the upper u level with the corresponding Einstein coefficients. 

where we assume no collisional (de-)excitation. The expression (pi^uiv) describes the fine profile of the transition of photons of 

poo 

frequency v and is normalized to one, that is, cf)(v) dv — 1. 

Similar to Eq. (B.21) we can write the extinction coefficient, which describes the radiative excitation from the lower to the upper 
level 

xTuv-^^iBi,uhu(.y\ (B.22) 


An 


where B^u describes the Einstein B-coefficient for extinction. 


In addition to spontaneous emission and extinction we have to take the stimulated emission into account, which can be described 
by adding a negative opacity contribution to Eq. (B.22): 


h V 

^/.w.y “ “j (BI Bl^ii ll^ Bii i^ 

An 


(B.23) 


where B^j represents the Einstein B-coefficient for stimulated emission. So, for < n„ B^j we get laser (maser) emission. 
The different Einstein coefficients are related to each other by the Einstein relations: 


_ D 

^LiJ ~ 9 

flight 

Si Bf^u ~ Su Buj- 


(B.24) 


Eollowing Eq. (B.2), the differential optical depth Ty is defined as 


h V 


dTy ^ Kyds (ill Biu - n« Buj) (^;,„(v) | ds 


! A 

light 4 I Su 1 JL / \ 
Au,i h/-n« 0 ;,h(v) 


Stt 


Si 


ds, 


(B.25) 


where we used the Einstein relations Eqn. (B.24) in the last fine. By assuming LTE conditions and therefore Boltzmann population 
distribution 

n.. P.. , ^ ^ ^ fr., , ,, ^ 

(B.26) 


^ ^ g(-Ey~E,)lkBB„ _ 8^ ^-hVyilksTc, 


ni Si 


Si 


we can rewrite Eq. (B.25) by using Eq. (B.20) 


dzy - 


r ^2 


nv^ \ / 


ds. 


(B.21) 


Einally, we have to integrate along the fine of sight and obtain the optical depth Ty 


8 7ry2 


(B.28) 


where Nu - f ds describes the column density of a certain molecule in the upper state. 
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In order to express Eq. (B.28) in terms of the total column density Atot = again with the Boltzmann population 

distribution 

oo oo oo 

Mot = X Z “ = ”> Z “ ^ • Q (rex), (B.29) 

!=0 i'=0 i=0 

where we used the partition function Q{T^^), which is dehned as sum over all states 


2(7’ex) = Z'^'^ 




/=0 


(B.30) 


For our purpose we rewrite Eq. (B.29) in terms of Nu and A^tot^ 

N„ 


g. 


tot = — e-“' —“ • QC^ex). 


Inserting this expression into Eq. (B.28) gives 

^2 


_ light 
“ 8 ;rv 2 


Mot g. 


e(rex) 






cr . , „ p-Ei/kaTa 

light . 

A^tot 


_ ^-hv..,/ksT^^ 0 ;.„(v), 


87rv2"“’"'‘“^ e(7’ex) 

where we used Eq. (B.20) to achieve Eq. (8) for a single line of a certain component and molecule. 


(B.31) 


(B.32) 
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